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Abstract 

This work aims at generation of optimum camber shapes for a wing for transport aircraft application. 

Induced drag reduction is the main consideration. The lift is increased by a factor and objective 
function for minimum induced drag is formed that is differentiable with respect to circulation that is 
aimed to generate prescribed lift. 
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Introduction 

Moderately swept back wing is considered for the application in favor of medium aspect ratio and low induced drag. 
Design and analysis is done for Mach number of 0.7 and 4° angle of attack. Several parameterization techniques are 
seen in literature that have been developed and applied over the years for aerofoil design; such as polynomial 
representation and non-uniform rational B-splines (NURBS). These offer precision in representing and manipulating 
analytical curves [1, 2]. In our work, forces acting on the wing are determined through panel numeric. Downwash 
influence coefficients are determined within the framework of linearised potential flow theory. Principles of calculus 
of variations are applied to form objective function [3, 4]. 

Constraint optimization technique is applied to the definition of camber lines for minimum induced drag for the 
specified value of lift. Program is developed in Fortran language. Gfortran complier is used that is available in Fedora 
operating system of VMware which has Linux based platform. Following command complies and creates a binary 
ready for executions. 

gfortran -o program name program name.f 

Following command executes the program file 
./program namecinput file name 


Following command opens the program in vi editor 
vi program.f 
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Following commands have been used in the gfortran 
complier. 

• pwd present working directory 

• cd change directory 

• mk dir create directory 

• vi program.f opens the program in vi editor 

• cpprog.f/directory name copy program 

• :q quit editor mode 

• se nu gives no sequence for entire file 

• rm program name removes a program 

• gfortran - o program name program name.f 
compiles and creates a binary ready for execution 

• ./program name executes the program file 

Encouraging results are obtained with this application. 


Mathematical Modeling 

Aerofoil shape parameterization can be obtained 
through definition of curves. NURBS is such a technique 
that utilizes a knot vector as defined below: 


C(f)=2 





In this formula ^ is the curve parameter, ranging from 0 
(the start of the curve) to 1 (end of the curve), k is the 
number of control points, E i n is the i-th basic function of 
order n, Wj is the weight associated with the i-th control 
point, and Pj=[Xj,yj] is the control point. 

In the numerical optimization, the change of aerofoil 
curve shape is achieved by changing the coordinates of 
NURBS control points. Control point in this approach 
cannot move along a direction more than a defined 
length with prescribed limits. Limits are to be imposed 
on one of the NURBS control points. Shape-change 
framework with Radial basis function (RBF) is another 
approach for aerofoil optimization. Using the RBF 
method, the aerofoil doesn't need to be parameterized, 
and is instead expressed as a cloud of points of arbitrary 
order and spatial resolution. There are a large number 
of basic functions to choose from. A radial basis function 
operates on radius between points, and returns a scalar 


value. The returned value will vary between 1.0 when 
the distance is 0, and 0 when the distance is equal to 
support radius. This support radius is chosen by the 
user, and represents the radius of influence of one point 
on other points. Camber changes occur through the 
camber cloud, which causes a change in the shape 
change cloud, which then in turn changes the external 
shape to reflect the change in camber. This essentially 
becomes trial and error technique. 

In the approach made herein, wing is represented by a 
large number of constant pressure panels to model 
circulation. These panels are used for estimation of 
pressure difference coefficient, lift and induced drag. 
Only half of the wing is considered. Other half of wing is 
imaged. The program uses vortex panel method to 
calculate downwash by satisfying tangential flow 
conditions to determine singularity strength of each of 
panels; and thus lift, induced drag and moments for the 
given planform are determined [3, 4]. Control point for 
determining downwash is taken at 0.95% of local panel 
chord. The imaging of half of wing is done though the 
following expressions: 

y = y k ~y P ( 2 ) 

yiimage) = -y k - y p 

The subscript p refers to panel and subscript k refers to 
control point location. Panel circulation (y) is 

determined through tangential flow boundary condition 
and pressure difference coefficient is given by 
ACp = 2y/U . Lift is determined through integration of 

A Cp over the chord and span. The program is 

developed in Fortran. It generates the output matrix of 
pressure difference coefficients from where lift, drag 
and moments are determined and then a optimization 
constraint of lift is introduced i.e., lift before and after 
the optimization remains same. Due to the effort of 
optimization, the given wing camber gets modified and 
hence a new pressure difference distribution results. 
Axis system and related moments are shown in Fig. 1. 
Half of wing is divided into nine chordwise and nine 
spanwise panels, making a total of 81 numbers of 
panels. 
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Figure l.Axis System and Related Moments 


Objective function (F) for the drag (D) minima and 
specified value of lift is written in the Lagrange form 

using Lagrange multipliers (\ ): 


F = D + A 0 (L - L) 


( 3 ) 


Where L is the specified value of lift which is 
expressed in the following way: 


L - Lx LIF 


( 4 ) 


Here L is the lift of a given profile that is to be increased 
by a factor LIF. 

Lift is expressed in terms of circulation as below for N 
number of singularities. 


L = pu[\y x +--- + A N r N ] 


Here. 7,. y N are circulation strength of N number of 

panels and A Av are re l ate d panel areas. Lift effects 

are simulated by a spread of vortex sheet (circulation 
strength). The property of the vortex sheet is that the 
component of flow velocity tangential to sheet 
experiences a discontinuity change across the sheet that 
is given by: 

7 = u i ~u 2 

Where u l and u 2 are tangential velocities just above 
and below the sheet respectively. 

Matrix of optimization Eq. (6) below is obtained by 
differentiation of objective function (F) w.r.t circulation 

y 


( 5 ) 
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Where XC x XC N are distances of panel control 


Trailing edge sweep = 0 .16 radians 


points from wing apex; and a i,j are influence 
coefficient towards downwash velocities. 

Solution of this matrix for prescribed lift L = LxLIF 
values results in optimal circulation from where 
pressure difference coefficient, drag, and wing root 
bending moment etc. are determined [5]. 

Results and Discussions 

Following wing is taken as candidate for study. Flow 
Mach number is taken as 0.7 and angle of attack is taken 
as 4°. Wing has practical application for high speed 
transport aircraft and as well as for UAVs. The value of 
LIF is taken as 1.2 for increasing the lift. 

Root chord = 2.9m 
Tip chord = 1.2909m 
Semispan=7.99m 

Leading edge sweep = 0.349 radians 


Flat surface wing shape is first considered and is 
subjected to optimization towards minimum induced 
drag. The lift is then incremented by a factor and again 
the optimal profile for minimum drag is worked out. 
Table 1 below shows the aerodynamic data comparison 
for the case of lift incremented wing with the base line 
wing profile. 

Table 1 shows that there is substantial potential for drag 
reduction through optimization. Figure 2 shows the 
values of chordwise pressure difference coefficient 
before and after optimization at root station for Case-A. 
There is shift in local lift towards the mid of chord. 
Similar feature is seen towards the tip of wing (Fig. 3). 
This is the reason for the increase in pitching moment 
ratio value. The increase in pitching moment is more in 
Case-B, since there is increase in local lift towards 
trailing edge. There is also shift in local loading towards 
tip as seen by the values of ratio of wing root bending 
moment before and after optimization (Case-A, value of 
1.0226 and Case-B, value of 0.8522). 


Table 1.Comparison of Aerodynamic Coefficient 


Conditions 

C L 

Before 

optimization C D 

After optimization 
Cd 

% Drag 
reduction 

M x /M x 

MJM y 

Without 

incremental 

value of lift 
(Case-A) 

0.381 

0.0266 

0.0093 

65% 

1.0226 

0.817 

After 

incremental 
value of lift by a 
factor of 1.2 
(Case-B) 

0.457 

0.0319 

0.0134 

58% 

0.8522 

0.6813 


Pressure difference coefficient for the Case-B is shown 
in Fig. 4 for two stations. There is similarity of shift of 
local lift with the Case-A. In both the cases, there is shift 


in location of adverse pressure gradient. The shift is 
towards the rear, and this is favorable from the point of 
delay in flow separation. 
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Figure 2.Pressure Difference Coefficient at Root Station without Lift Increment (Case-A) 
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Figure 3.Pressure Difference Coefficient at Tip Station without Lift Increment (Case-A) 
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Figure 4.Pressure Difference Coefficient with Lift Increment (Case-B) 
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Root station 



Tip station 



k/c 

Figure 5.Camber in Radians after Optimization before Lift Increment (Case-B) 


Root station 


Tip station 



Figure 6.Camber in Radians after Optimization after Lift Increment (Case-B) 


Resulting optimal camber is shown in Fig. 5 and 6 in the of spanwise twist is lesser in Case-B. This is visible by the 

chordwise manner (note that the initial camber in both lesser displacement of camber values from root to tip. 

the cases is nil, i.e. a flat surface case). The requirement 



Root station, Chordwise 



Tip station, Chordwise 



Twist .Spanwise 
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The optimal warp of Case-B is separated for camber and 
twist by a using a subroutine given at Appendix 'A'. 
Resulting values of twist and camber are shown in Fig. 7. 
The resulting aerofoil camber is seen to be same at root 
and tip. This is very advantageous from manufacturing 
ease point of view. The twist is also nearly linear. 


Excellent features of design are obtained by ths 
technique followed herein. 

Values of % camber in chordwise manner are worked 
out and shown in Fig. 8. A single aerofoil results at root 
as well at tip across the span after separation of twist. 


* root as well at tip after separating twist 
root before separation from twist 
- tip before separation from twist 
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Figure 8.Values of % Camber in Chordwise Manner 
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Appendix 'A' 

C WING WARP IS BROKEN INTO TWIST AND AEROFOIL CAMBER 

C N=NUMBER OF PANELS, XP= PANEL CORNERS , RT=TWIST 

C AZ= GIVEN WARP TO BE BROKEN INTO AEROFOIL CAMBER + TWIST 

SUBROUTINE TWISTR(XP,ZZ,RT,AZ,CAMBER,N) 

DIMENSION XP(10,10),ZZ(10,10),RT(10),AZ(9,9) 

DIMENSION XZZ(10,10),CAMBER(9,9) 

DATA IR,IW/5,16/ 

N1=N+1 
DO 100J=1,N 
DO 150 I1=1,N1-1 
l=N-ll+l 

DIS=XP(I+1,J)-XP(I,J) 

150 ZZ(I,J)=ZZ(I+1,J)-DIS*AZ(I,J) 

100 CONTINUE 

WRITE(IW,200) 

200 FORMAT(/2X,'VALS OF THE VER-ORDI ROWS READ CHORDWISE’) 
WRITE(IW,225)((ZZ(I,J),I=1,N),J=1,N) 

DO 250J=1,N 

250 RT(J)=-(ZZ(1,J)-ZZ(N+1,J))/(XP(N+1,J)-XP(1,J)) 

WRITE(IW,275) 

275 FORMAT(/2X,'VALS OF THE SPANWISE TWIST DISTRIBUTION') 

WRITE(IW,225)(RT(J),J=1,N) 

DO 300 J=1,N 
DO 350 1=1,N 

350 CAMBER(I,J)=AZ(I,J)-RT(J) 

300 CONTINUE 
WRITE(IW,375) 

375 FORMAT(/2X,'VALS OF THE AEROFOILS CAMBER DIST') 
WRITE(IW,225)((CAMBER(I,J),I=1,N),J=1,N) 

DO 660 J=1,N 
XZZ(l,J)=0. 

DO 677 1=2,N 

677 XZZ(I,J)=XZZ(I-1,J)+CAMBER(I,J)*(XP(I+1,J)-XP(I,J)) 

660 CONTINUE 
DO 899 1=1,N 
DO 899 J=1,N 

899 XZZ(I,J)=XZZ(I,J)/(XP(N+1,J)-XP(1,J)) 

WRITE(IW,775) 

775 FORMAT(/2X,'VALS OF AEROFOIL NON-DIM 

1VER-ORDI ROWS READ CHORDWSE') 
WRITE(IW,225)((XZZ(I,J),I=1,N),J=1,N) 

225 FORMAT(9F8.4) 

STOP 

END 
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Appendix 'B' 

Nomenclature 
A = Panel Area 


a i j = Influence coefficients (Influence of j th panel on i th control point) 
b = Span 

c = Local chord length 

C D = Drag coefficient 

C L = Lift coefficient 

ACp = Pressure difference coefficient 

D = Induced drag 

L = Lift 

LIF = Lift Increment Factor 

M = Mach number of freestream 

N = Total number of panels on half of wing 


M x = Wing root bending moment before optimization 

* 

M x = Wing root bending moment after optimization 

M = Pitching moment before optimization 

* 

M = Pitching moment after optimization 


x, y, z = Cartesian coordinates in geometric plane ; x is chordwise, y is spanwise and 
z is vertical 

XC = x-Distance of panel control point from the apex of wing 
dz/dx = Aerofoil camber slopes (ZX) 

U = Freestream velocity 
a = Angle of attack 



y = Circulation 
p = Density 

Suffix 


i = i th Panel control point 


j = j t h pane| 
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